//
//  BMvanilla.cpp -- this is a plain implementation of the Bratley-McKay algorithm  
//  This was to test that my description of the BM algorithm to Raghavendra 
//   was correct. It turns out that this is a lot cleaner.
//  
//
//  Created by Suresh Govindarajan on 27/07/12.
//  Copyright (c) 2012 suresh.govindarajan@gmail.com. All rights reserved.
//
#include <cstdio>
#include <iostream>
#include <cstring>

using namespace std;

#define FOR(i, a, b) for(int i=int(a); i<=int(b); i++)
#define REP(i, n) for(int i=0; i<int(n); i++)
#define RFOR(i, a, b) for(int i=int(a); i>=int(b); i--)

const int NMAX=100,MAX=NMAX;
typedef  unsigned long ll;
//
// A node is stored as a single long integer with SHIFT bits per coordinate. 
// Thus if dim=4, 40 bits are used by the node -- we can go up to 6 dimensions
// with the current SHIFT=10.
//
// A Partition is a collection of nodes
//
typedef ll Partition[NMAX];

//
// Choose dim=2,3,4,5,6 
// 2=ordinary partitions; 3=plane partitions; 4=solid partitions
//
const int dim=4;
const ll  SHIFT=10;
const ll  ONES = (1<<SHIFT)-1;
//const ll ELBIT = -1-ONES;
//const ll ELBIT=(ONES<<(SHIFT))+(ONES<<(2*SHIFT)); // erase last coord

/************** GENERATING TOPOLOGICAL SEQUENCES *****************/
ll d[MAX+1];		// d[i] = no of topological sequences of index i
//int cc[MAX+1][MAX+1];	// no. of columns in given plane,row
//int pp, n;		// no. of planes and n
int n;
//const ll dxx=(ll(1),ll(1)<< (SHIFT),ll(1)<< (2*SHIFT));
Partition current;
Partition ipossible;
const ll dy[6]={ll(1),(ll(1)<< (SHIFT)),(ll(1)<< (2*SHIFT)),(ll(1)<< (3*SHIFT)),(ll(1)<< (4*SHIFT)),(ll(1)<< (5*SHIFT))};



//typedef pair<int,int> pii;
//typedef pair<pii, int> piii;

// increment the i-th coordinate by one
inline ll dx(const int &i) {
	return (ll(1) << (i*SHIFT));
}

// get the ith coordinate
inline ll coord(const ll &x, const int &i) {
	return (((x>> i*SHIFT) & ONES));
}

void printNode(ll node) {
	RFOR(i,dim-1,0) printf("%ld ",coord(node,i));
    printf("\n");
}

void printPartition(Partition node,int size) {
	printf("( ");
    FOR(j,1,size){
        RFOR(i,dim-1,0) printf("%ld ",coord(node[j],i));
        if (j!=size) printf(", ");}
    printf(")\n");
}

//checks if the presented node is in partition
bool inpartition(Partition partition, ll node, int size){
    int i;
    FOR(i,1,size) if(partition[i]==node) return true;
    return false;
}


// Checks if a given node can be added in a given partition -- it looks to see if the supporting nodes are present.

bool ispossible(Partition partition, ll node, int size){
    int j;
           FOR(j,0,dim-1) if ((coord(node,j)!=0) && (inpartition(partition, (node-dy[j]),size)==false)) return false;
    return true;
}


// a[i] = number of 1-d partitions of sum i //
ll a[MAX+1];
void calculate_1d(int n) {
  for(int i=0; i<=n; i++) a[i]=1;
  for(int k=2; k<=n; k++)
    for(int i=k; i<=n; i++)
      a[i] += a[i-k];
}


// The main recursive routine for the code

void part(int size, int from, int to){
    int i,j,k,bean=0;
    ll tmp;
    FOR(i,from,to){
        current[size+1]=ipossible[i];
        d[size+1]++;
        if (size==n-1) {
            //printPartition(current,size+1);
            continue;
        }
        
        REP(j,dim){
            tmp=ipossible[i]+dy[j];
            if(ispossible(current,tmp,size+1)==true){
                  bean++;
                    ipossible[to+bean]=tmp;
                }
        }
                //printf("bean=%d ",bean);
        //printPartition(ipossible,to+bean);
        part(size+1,i+1,to+bean);
                            bean=0;
            }
}


void generate_topological_sequences(int nn) {
    d[0] = 1;d[1]=1;
    FOR(i,2,NMAX) d[i]=0;
    //FOR(i,1,dim*NMAX) ipossible[i]=0;
    FOR(i,1,NMAX) current[i]=0;
    //current[3]=ll(1);current[2]=ll(1)<<SHIFT;//current[4]=ll(3);
    //current[2]=ll(1)<<(2*SHIFT);current[3]=ll(2)<<(2*SHIFT);
    ipossible[1]=ll(1);
    FOR (i,2,dim) ipossible[i]=(ll(1)<<((i-1)*SHIFT));
    //ipossible[2]=(ll(1)<<SHIFT);ipossible[3]=(ll(1)<<(2*SHIFT));//ipossible[4]=ll(1)+(ll(1)<<SHIFT);
      n = nn;
        part(1,1,dim);
}

//////////////////////////////////////////////

int main() {
    //printf("Enter N: ");
    int N=15; //scanf("%d", &N);
    //generate_topological_sequences(N);
    // Initialisation
    d[0] = 1;d[1]=1;
    FOR(i,2,NMAX) d[i]=0;
    FOR(i,1,NMAX) current[i]=0;
    //current[3]=ll(1);current[2]=ll(1)<<SHIFT;//current[4]=ll(3);
    //current[2]=ll(1)<<(2*SHIFT);current[3]=ll(2)<<(2*SHIFT);
    ipossible[1]=ll(1);
    FOR (i,2,dim) ipossible[i]=(ll(1)<<((i-1)*SHIFT));
    n = N;
    // The most important call -- the current partition has one node and ipossible contains dim nodes to start with.
    part(1,1,dim);
    printf("d[i] = No. of %d-dimensional partitions of i\n",dim-1);
    for(int i=0; i<=N; i++) printf("d[%d] = %lu\n", i, d[i]);
}


